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ABSTRACT 


The investigation of the vibrations of a turbine blade 
has been a fascinating and challenging problem to the design 
engineers and the research workers and they have 

i) idealized it as a beam, a plate and a shell, 

ii) incorporated one after another the additional effects of 
asymmetry of the cross-section, taper, pre-twist, stagger 
camber, packetting, rotary inertia, shear deformation 
etc* and 

iii) utilized different appro ximate/nisnerical techniques Of 
solution like Rayleigh-Ritz, Ritz-Galeikin, Holzer, 

A 

Collocation, FDM, FEM etc* 

t 

Here, the free vibration of a stationary blade has been 
dealt with modelling it as a beam and as a plate* The finite 
element technique using Galerkin* s method has been employed to 
reduce the system of partial differential equations to a matrix 
eigen-value problem* The eigen-values and the eigen-vectors 
were computed using a NAG subroutine F02AEF which makes use of 
the Householder* s decomposition and the QL algorithm* As the 
first eigen-value, thus computed, exhibited oscillating trend 
rather than converging to its lower bound with the increase in 
the number of elements, in lieu of the eigen-value problem, an 
Inverse eigen-value problem was solved* 



To determine whether a blade of asymmetric cross-section 
can be treated as a beam, one needs to compute the chord length 
and the maximum thickness of the cross-section* This has been 
done invoking Zukovsky* s transformation* While treating the 
blade as a beam, the asymmetry of the cross-section and the 
pre-twist accompanied by the fibre-bending during torsion were 
accounted for* The frequencies were computed for three differ- 
ent data of the blade* The effect of pretwist on the coupled 
bending-bending modes and the torsional modes for a beam of a 
rectangular cross-section was studied. When treated as a plate, 
the blade was assuned to be so thin that the effect of the 
rotary inertia and the shear deformation could be neglected* 

In computation of the element stiffness matrix and the element 
mass matrix, integration was carried out by a Gauss quadrature 
formula. The frequencies were computed for four different 
aspect ratios of the blade* All the results obtained were 
compared with those obtained by the respective authors* 



CHAPTER 1 


INTRODUCTION 

The turbomachlnes are widely used in the steam tur- 
bines, gas turbines, compressor and gas turbines of the air- 
craft propulsion systems, marine drives etc* The steam turb- 
ines comprise of three stages called as high pressure (HP), 
intermediate pressure (IP) and low pressure (LP)j each stage 
consists of the moving blades fixed to the rotor and the stat- 
ionary blades fixed to the casing* The turbine blades are the 
soul of the whole turbomachine* The minimum satisfactory height 
of impulse blading is from 0*1 to 2*5 cm in small turbines and 
from 4*0 to 5.0 cm in large turbines [23] • The values of exit 
angle for the turbine blades in vogue vary from 150° to 165° 
at HP and IP stages and from 140° to 150® at the LP stage* 
sometimes reaching 130° to 140® in the final stage of large 
turbines where maximum flow is needed* The values of the 
chord length of the bladings of the impulse turbine range from 
4*0 to 5*0 cm at the HP end to 15 cm or more at the LP end* 

The HP turbine blades are short, and tapered, have nearly equal 
entrance and exit angles thus requiring small angle of the 
pre-twist — say of the order of 2®, have stagger angle between 
25® to 35®, and are banded at a definite pitch. The LP turbine 
blades are long with considerable pre- twist. In general, a 
turbomachine blade is a tapered, pre- twisted, staggered and 
cambered cantilever of asymmetric airofoil cross-section* 
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Because of asymmetry, the centroid and the flexural centre Of 
the aerofoil section do not coincide# 

The working fluid, throttled through a noirle, flows 
over a blade and in due course exerts a tangential force on 
it, thereby causing it to rotate about the turbine axis and to 
undergo bending and torsional deformation* The process of 
power generation in the turbine leads to the excitation of the 
blades because of the flow pattern of the working fluid. 

Cyclic stresses are set in the blades due to forced vibrations. 
They are detrimental to the blades as they are likely to cause 
fatigue failure* Moreover, if the forcing frequency is close 
to the natural frequency of the blade, it leads to the pheno- 
menon of resonance thus further enhancing the amplitude of 
cyclic stresses and hence increasing the probability of fatigue 
failure of the blade. Therefore, it is imperative to know the 
natural frequencies of the turbine blade* 

The vibration analysis of the turbine blades has been 
evolved from using the simple beam theory to that of the 
sophisticated finite element techniques. When idealised as a 
simple cantilever beam of symmetric cross-section (rectangular), 
the analysis of vibrations of the blade is quite simple* The 
analysis of the blades with aerofoil section has been arrived 
at, Incorporating successively the effects of taper, asymmetry, 
pre-twist, fibre-bending, Coriolis force, centrifugal force, 
stagger angle, rotary inertia and shear deformation# Additional 
effects due to factors like residual stresses, relief of 
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stresses of cambered section during bending, the penetration 
of the stress pattern into the root, the angle at the root and 
tip to the axis of the blade, non-coincident flexural and tor- 
sional centres etc* have been studied by various researchers • 
Modelling the blade as a cantilever beam with asymmetric 
aerofoil section, pre-twist and taper can give good results 
when the length of the blade is much greater than the other 
dimensions* However, a short flat blade has to be idealized 
a* a thin plate and a short twisted blade a thin shell* The 
effects of rotary inertia and shear deformation have also been 
incorporated and studied by many researchers and scientists 

ua.os# m * 

When idealized as a uniform straight beam of asymmetric 
atrofoil section, the blade exhibits the coupled vibrations of 
flap-wise bending, chord-wise bending and torsion* This is 
owing to the asymmetry of the section by virtue of which the 
eenttoid and the flexural centre are non-coincident whereas 
when considered as a beam of rectangular cross-section (the 
centroid and the flexural centre being coincident) coupled 
bending-bending vibrations occur which are attributed to the 
pre->»twlst* Because of pre-twist the product of Inertia comes 
into picture which couples the flap-wise bending and chord- 
wise bending | but the torsion remains uncoupled* The only 
factor to couple torsional modes to the bending modes is 
distinct location of centroid and flexural centre* If the 
blade is assumed to behave as a flat rectangular cantilever 
plate, no coupling occurs between bending and torsional modes. 
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as the displacements along the axis and the chord are directly 
proportional to the slopes along the axis and the chord res* 
p actively* Nevertheless, in the case of the plate also, 
coupling occurs only if it is pre-twisted in which the two 
dimensional element has to be treated as a shell* In all the 
cases discussed above, it is difficult to obtain the closed 
form solution* Hence one has to resort to the approximate/ 
numerical methods like Rayleigh-Ritz, Holzer, Collocation, 
Galerkin, FDM, FEM etc. 

BRIEF REVIEW OF THE AVAILABLE LITERATURE ! 

Houbolt and Brooks Qf] developed the differential 
equations of motion for the lateral and torsional deformations 
of twisted rotating beams without the loss of generality as 
regards the non-coincidence of the torsional, flexural and 
centroidal axes for the application to helicopter rotor and 
propeller blades* Carnegie [[£] derived the governing equations 
of coupled bending-bending* torsion motion of thin pre-twisted 
cantilever blade using variational method* He incorporated 
the additional effect due to pre-twist while studying the 
vibrations of a straight, pre-twisted uniform blade of rect- 
angular cross-section as a special case* Paranjpe l. 5] studied 
the self-sustained oscillations of the turbomachine blading 
which were likely to occur particularly when light sections 
were subjected to heavy aerodynamic loading* Carnegie Q4H ~4 

derived the general equations of motion of the rotating canti- 
lever blade applying the variational method. Kulkami et al {jsH 
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applied both the Reisner method and the total potential 
energy approaches to deduce the equations of the coupled 
bending-bending vibrations of the p re- twisted uniform rect- 
angular cantilever blade accounting for the shear-deformation 
and the rotary inertia effect* They found out the natural 
frequencies by Ritr process and compared these with the 
experimental results for a specific data* 

Duggan and Slyper Q5] applied a modified central 
finite difference method to the equations of motion to compute 
the natural torsional frequencies of vibration of the pre- 
twisted cantilever beams of non-uniform symmetric cross-section 
considering the effect of fibre-bending* Belgatankar et al l6j 
obtained the natural frequencies of the torsional vibration 
of a cantilever beam of rectangular cross-section with uni- 
form taper by collocation method and compared these with those 
determined by Holzer* s method [_<£]+ Carnegie and Rao H7j app- 
lied Galerkin* s procedure to determine the natural frequencies 
of an untwisted cantilever blade with an asymmetric aerofoil 
cross-section executing coupled bending-bending- torsion 
vibrations, neglecting the effect of the fibre-bending* Rao 
f8] computed the coupled bending-bending frequencies of a 
pre- twisted cantilever blade of the uniform rectangular cross- 
section* He determined the coupled bending-bending- 
torsion frequencies of vibrations of a uniform untwisted 
cantilever beam of an asymmetric aerofoil section using the 
collocation method and ignoring the effect of the fibre- 
bending* He D® computed the natural frequencies of the 
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torsional vibrations of a pre-tvdsted cantilever beam with 
the uniform rectangular cross-section. 

Carnegie and Rao fif] found out the natural frequencies 
of the torsional vibrations of the tapered y pre-tvdsted 
cantilever beams of the rectangular cross-section. Gupta 
and Rao S.S. P13J used finite element method (FEM) to study 
the vibration of a single, tapered, pre-tvdsted blade of an 
aerofoil cross-section* Ba j a j and Dhoopar Q 14 j studied the 
vibrations of the packetted, tapered blade of the rectangular 
cross-section with the negligible pre-tvd st applying FEM. 

Kara and Subrahmanyam Cl Cl solved the problem of calculating 
the natural frequencies and the mode shapes of the rotating 
blades by a finite difference method (FDM). Saurabh Kumar ClC 
tackled the- problem of the coupled bendin g-bendi ng-to rsi onal 
vibrations of a single cantilever blade with the effects of 
the pre-tvd st, rotation, hub- radius, stagger angle and the 
asymmetry of the cross-section by FEM using the Galerkin* s 
technique* Rao [JSf] compiled the results pertaining to the 
effects of the speed of rotation, asymmetry, stagger angle 
and hub radius on the natural frequencies of the vibrations 
of the pre- twisted cantilever bladings used for the turbo- 
machines. 

In their state of art paper, Leissa, Kielb and Macbain 
[XT] compiled the results of vibration of the twisted canti- 
lever plate obtained by a joint govemment/industry/universlty 
research study. The numerical results were obtained for 
twenty different geometries using nineteen different methods — * 
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15 using FEM f 2 using shell theory, 2 using beam theory — in 
order to have comprehensive, definitive set of results. The 
graph of the non-dimensionalised frequency versus the pre-twist 
exhibited different trends in some cases, so much so, that the 
predicted dependency of the frequency on the twist-angle for 
the twisted cantilever plates spanned the spectrum from incre- 
asing significantly to decreasing significantly. This, indeed, 
made the problem quite challenging and fascinating. 

SCOPE OF THE PRESENT WORK : 

The present work aims at studying the free vibration of 
the turbine blade treating it as a) a beam and b) a plate by 
finite element method (FBI). It is observed that the effect of 
rotation of the blade is negligible, hence we study the vibra- 
tions of non— rotating blades. To limit the study to a manage- 
able size, the effects of stagger, camber, taper, centrifugal 
force, Coriolis force, rotary inertia, shear deformation, packet- 
ting, hub-radius etc. are neglected. 

As reviewed above, vibrations of the blade treating it 
as a beam have been analysed by various authors using classical 
approximate methods. The objective of this work is to re-analyse 
some of these problems using FEM thereby establishing that FEM 
can be used sufficiently accurately to obtain the natural frequ- 
encies of straight or pre- twisted long blades by treating it as 
a beam. FEM has been used by some authors earlier also but not 
so extensively. It is also established that a short straight 
blade has to be analysed as a plate and here also FEM can be 
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used with reasonable accuracy. A short twisted blade has to be 
modelled as a shell. Due to lack of time, this study could not 
be undertaken. 

Chapter 2 deals vdth the problem of vibration of a blade 
considering it to be: 

1) straight, untwisted beam of asymmetric aerofoil section in 
coupled bending-bending-torsional modes caused due to the 
asymmetry of the section. The data used in the reference 
[l2] to solve the problem is used to compute numerical 
values of the natural frequencies and to compare the results 
with those obtained in the same paper. Coupled modes of a 
straight cantilever blade are studied assuming that the 
principal axes of the aerofoil section are closely located 
to the reference axes, say of the order of 2°. In other 
words, the product of inertia with respect to the reference 
axes is assumed to be negligible. The effect of fibre- 
bending is ignored in this case. 

2) Straight beam of rectangular cross-section — untwisted as 
well as pre- twisted. For untwisted blade, flapwise bending, 
chordwise bending and torsion are all uncoupled. In the 
case of pre-twisted beam, flapwise and chordwise bendings 
get coupled due to pre- twist but torsion remains uncoupled 
as centroid and shear-centre of the section are the same. 
Using the data of the references [18] and [8] and incorpor- 
ating the effects of pre- twist and fibre-bending, the 
results are obtained by FEM and compared with those given 
in the same paper. Using the data of the reference [17] 
and ignoring the effect of fibre-bending, the natural 
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frequencies are computed and the results are compared with 
those compiled in this paper. 

Chapter 3 consists of modelling the blade as a simple canti 
lever uniform thin plate of rectangular cross-section. The flap- 
wise bending, chordwise bending and torsion are not coupled becaus< 
displacements along axis and chord are proportional to the slopes 
along the axis and the chord respectively. In this case, the data 
from the reference [17] is used and the natural frequencies obta- 
ined are compared with those given in the same paper. 

In each case, as the system of partial differential equa- 
tions is coupled and its closed form solution is difficult to 
obtained mathematically, the system of equations is solved by Finit 
Element Method. This comprises of discretizing the blade into a 
finite number of elements, approximating the deflections over each 
element by suitable shape functions, computing element stiffness 
and mass matrices thereby converting the system of the partial 
differential equations to a matrix eigen-value problem, assembling 
the element stiffness and mass matrices, applying the essential bou 
dary conditions and solving the eigen— value problem. Finite elemen 
method, as described above, is employed because results obtained by 
this method converge to the exact solution with the refinement of 
the mesh of elements. Two noded beam elements and four noded plate 
elements have been employed in Chapters 2 and 3 respectively. To 
solve the eigen value problem, a NAG subroutine titled F02AEF [30] 
is made use of. This subroutine uses Householder's decomposition 
and QL algorithm to compute eigen-values and eigen vectors. The 
graphs of Chapter 2 have been plotted using GPGS graphic subroutine. 
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CHAPTER 2 

VIBRATIONS OF BLADE AS A CANTILEVER BEAM 

2.1 BLADE AS A BEAM ? 

When the length of the blade is large as compared to 
its other dimensions, it can be treated as a beam. This 
chapter deals with the problem of free vibration of the stat- 
ionary turbine blade idealized as 

a# untwisted straight beam of asymmetric aerofoil cross- 
section and, 

b. untwisted as well as p re- twisted straight cantilever 
beams of rectangular cross-section. 

The effects of stagger angle, taper, rotary inertia and shear 
deformation* are neglected. 

In the case of an asymmetric aerofoil-section, its 
centroid and the centre of flexure are distinct. This give* 
rise to coupled bending-bending- torsional vibration even when 
the pre-twist of the beam is sero. The turbine blading ha# 
asymmetric aerofoil section with one principal area moment 
of inertia much greater than the other. 

For an untwisted uniform beam of rectangular ero##- 
section, it is observed that the flapwise bending, the chord- 
wise bending and the torsional modes of vibration are all 
uncoupled. But, in the case of a pre- twisted beam with rect- 
angular cross-section, bending-bending coupled motion occur# 
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due to p re- twist and torsion is still uncoupled as the cent- 
roid and the flexural centre are coincident* 

2*2 DERIVATION ? 

Consider a pre-tvdsted straight cantilever beam of an 
asymmetric cross-section as shown in the Figure 2*1. Let 
u t v and © be chordwise displacement, flap-wise displacement 
and angular displacement about the axis of torsion of the beam 
respectively* The strain energy \jf\ of this beam is given by 
the expression: 

V “ 2 o ; f E C I yy ( 0> + 2I xy ( f^ )( ^ ) + I xx ( ^ ) > 

G< ff )Z - C 2 ( f^> 2 } * 

where I » I and I are the area moments and product of 
yy ** xy 

inertia of the cross-section of the beam about the reference 
axes* The term containing I is because of the pre-twist r 

Xy 

and it causes coupling of flapwise and chordwise bending. 

Let I an <* be the principal area moments of inertia of 
the cross-section. The area moments and product of inertia 
depend on the angle between the reference axes and the prin- 
cipal axes of the section and their transformation is given by 

I ** I** 1 + I . A. (see reference L& ) 

« 9V « m *6 


where 
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I yy I xy 


I I 
xy xx 


A u 

*b 


cos 


sin 


2Tz . 2 Tj 

t 8111 nr 


2^2 


--cos 


ZlM 


av 



and I 


am 


f.2 2 : , . x a i, 
2 


(2*1.b) 


C is the torsional stiffness of the cross-section and is given 
by 

2 


r * 

c * C 1 + C 3^ 


( 2. 1 . C ) 


The term containing C 2 in equation (2*1 .a} is the contribution 
of fibre-bending during torsion* 

and C 3 are the functions of material properties 
and geometry of the beam. The terms corresponding to these 
constants are due to the total torque on the cross-section 
that comprises of 

i) St* Venant torsional moments This is because of the 
twisting as well as warping of the cross-section, 
ii) torque due to pre-twist: The pre-twist augments the 
warping of the cross-section thereby increasing the 
restoring torque and 

ill) torque due to fibre-bending: The tendency of the 

beam to reduce warping thus trying to retain the 
planeness of the cross-section during the deformation 
results in the additional torque which is the conse- 
quence of the bending of each individual fibre* 
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where r and r are the co-ordinate distances between the 
x y 

centroid and the centre of torsion (which is assumed to be 

coincident with the flexural centre of the cross-section). 

The term containing r and r in the equation (2.2.a) is due 

* y 

to asymmetry of the cross-section which is responsible for the 
coupling of torsion and flapwise as well as chordwd.se bending. 
If the cross-section is symmetric, say rectangular, the cent- 
roid and flexural centre are not distinct and the torsional 
mode is decoupled from the bending modes# I c ^ Is the polar 
moment of inertia about the centre of flexure and is given by 

ht “ I cm +A(r x +r y ) (2 * 2>b) 

where I is polar moment of inertia about the centroid of 
cm 

the cross-section* 

Applying Hamilton’ s principle, 

t 2 

S / (T - v)dt «* 0 , (2. 3. a) 

*1 

One obtains the governing equation of coupled bending-bending- 
torsion vibrations as 
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[ A «*]' +| I* + C H ** 0 (2.3#b) 

in which dash {* ) and dot (*) represent partial derivatives 
with respect to spatial co-ordinate z and time t respectively* 
and 


T 

U = £u V 

9] , A * 

El 

y y 

EI 

xy 

0 

Eljfy 

EI xx 

0 

0 

o 

- C 2_. . 



«5 A 

0 

?Ar y 

B * diag \0 

R 

on 

Ok 

*oi 

» 

o 

0 

?A 

-?Ar x 



<?Ar y 

~ ?Ar x 

? X cf 


(2.3.c) 

It is obvious that tensor A causes bending-bending coupling 
whereas C causes bending-bending- torsion coupling* 

W 

At the fixed end, the deflections and the deflection 
derivatives are zero* Hence, the geometric boundary condi- 
tions are 

at z *= 0, u = 0 and u* = 0 (2«3*d) 

At free end, the shear forces, the bending moments and the 
twisting moment are all zero* Hence, the natural boundary 
conditions are 

at z •* 1 , A ju* * 0 

and ["A u«] * + B u 1 =0 


(2*3.e) 



First of all, consider an untwisted straight beam of 
asymmetric aerofoil cross-section in which case the effect 
of fibre bending is negligible i*e* # for X * 0, C 2 is small 
and hence neglected and 

A * diag[HI 22 EI^ Oj (2.3*f) 

The tensor C remains unaltered resulting in bending-bending- 
torsional coupled mode of vibration* 

As closed form solution of the system of equations 
(2#3*h) # (2.3*d) and (2.3*e) is difficult to obtain, to solve 
it, one has to apply one of the approximate numerical methods 
such as Gal exkin. Collocation, FDM, FEM, etc* Since the 
results obtained by FEM converge to the exact solution with 
the refinement of the mesh, FEM is utilised to solve this 
system of equations# 

For this case, since C 2 is small, in the variational 
functional of the equation (2.3*a), the highest derivatives 
of u and v are of the second order and that of 0 of the first 
order. Hence, the principal derivatives are (i) and u, 

(ii) and v and (ili) 9. The compatibility condition states 
that the approximating functions should be such that u, 

«• tl' 9 8houW bs <*»>«•.»»«. across th « loterelsmont boundary. 
The highest derivative occurring in the variation is of the 
order three for u and v and of the order for one 0 (as seen 
from the natural boundary conditions, equation (2.3.e)). 

Hence, the completeness condition implies that the shape 



t6 


function* for u and v should be cubic polynomial* and those 

for 9 should be linear polynomials. 

The beam is discretized into a finite number of two 

noded elements# A typical element is shown in the Figure 2*2# 

at*. 

There are five degrees of freedom at each node. Let u, , rrt 

3V 1 1 

v. » 5— » ©. be local degrees of freedom at one of the nodes 
1 au 0 av 0 

of the element and u 2 , Xz * v 2 , © 2 be tilose another 

node. One can then express u at any point P on the element as 


_ tT e i©t 
u « S U e 


(2.4*a) 


where 
(U®> #* 


au 

L u i 5z 


av 


v 


1 az 


e 


0 J 


2 32 


2 az 


'2 


and 




sSS 


P t P 2 0 ° 0 P 3 P 4 ° ° 0 

0 0 P t P 2 0 0 0 P 3 P 4 0 

0000 Q t 0000 Q 2 


( 2.4# b ) 


Here P^ s (i = t to 4) and Cl’ s (j = 1 and 2) are the func- 
tions of the local co-ordinate ?• 

For e element, the relationship between X and the 
global co-ordinate z is given by 

* a (* - z c ) 


where 



z 


e 


a * 
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is the half-length of the element. 



is the global co-ordinate of the mid-point and 
z g and are the global co-ordinates of the nodes. 

(2.5) 

The expressions for PM s and Q ^ 1 s which satisfy the 
convergence criteria are 

P 1 = 1 (2 - 31 + S 3 ) 

P 2 = | (1 - I - f + t 3 ) 

P 3 = | (2 + 3i - i 3 ) 

( 2 . 6 ) 

P 4 = f (- 1 - t + I 2 + i 3 ) 

Q, = 5 (1 - £) 

Q 2 = 5 (1 + S) 

Substituting u = S U e as an approximate solution 

in the equation (2.3.b), one gets the residue 

e e = {[4 i”] « + B y e - a 2 g s u e (2.7. a) 

B 

£ is a non-zero vector albeit it tends to zero as the size 
of the element is reduced. 



t8 


Applying Galericln* s method, on© gets 


‘e +1 


f 

z 


S T £® dz 


(2*7*b) 


© 


whence the system of partial differential equations simplifies 
to the following Eigen-value problem 


K e U e - A M e y e 


.e 


( 2*7 *c) 


where * and the expressions for the element stiffness 
matrix, the element mass matrix and the element stress-resul- 
tant vector are as follows* 


K 


e 


z ©+1 

; [ <i T f (A)(S)" - (sV (Bi(sy]dz (2.7,4) 


M e * 

as 


and 


R 


04-1 * 1 * 

2 # jk 

0 


F- S T (A U*)» + B u» + (S T ) , (A u«)} 


'e+t 


e 


(2*7*©) 


(2*7*f) 


To evaluate the element stiffness and mass matrices, 
one changes the variable from the global co-ordinate z to the 
local co-ordinate £ given by the equation (2*5)* Then 

dz ** a dx , 

3 | 1 3 ^ 

~ and ** ^*~2 


(2.7.g) 
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and the limits of the integration change from z e and z^ to 
-1 and 1 respectively, thus 

*e+1 1 

/ f(z)dz m a f f(z + a? )df (2*7*h) 

z a -1 

e 

After integration, one obtains element stiffness matrix and 
element mass matrix given by the equations (2*7*i) and (2*7*j ) 
respectively. These are symmetric po sitive definite matrices 
of the order 10 x 10. The element matrices, thus obtained, 
are assembled together to obtain the global stiffness matrix 
and the global mass matrix according to the standard procedure 
of FEM and each of them is singular. 

To solve this system of equations, one partitions the 
global matrices after applying the geometric boundary condi- 
tions 

dU JV 

at z*0 # Uj * 0, ^2 ** 0, V.j * 0, .g-j—" * 0, 0.J * 0 

(2.7.k) 

and hence one gets an eigen-value problem as 

KU-'AMU = 0 (2.7.1) 

33X mm SC *** *** 

Solution of this equation gives the eigen-values s in 
number equal to the order of K or M and the corresponding 
eigen-vectors. 




(2.&U 


($*L m Z) 



22 


2,3,! SOLUTION i 

To solve the equation (2.7,1) for the eigen-values 
and eigen-vectors, the NAG subroutine F02AEF D3QJ has been 
utilised. According to this ti for the system cf equations 
[a] { x } where [Aj must be a fjpmetric square matrix 

and [B] must be a positive definite symmetric^square matrix, 
all the eigen-values and eigen-vectors are computed using 
Ho us bolder* s decomposition and QL algorithm. Note that the 
eigen-values ^ represent the square of the natural frequen- 
cies o) 4 of vibration of the beam. As the continuous system is 
approximated by a set of fipite number of elsments, the infinite 
number of degrees of freedom corresponding to the continuous 
system is reduced to a finite number, thereby making the 
approximated model stiffer than the actual system, .Obviously, 
with the increase, in the number of elements, the number of 
degrees,of freedom increases, the overall stiffness of the 
approximated model decreases and consequently the eigen— values 

converge to their lower bounds. 

During the computation of the results, it was exper- 
ienced that as the number of elements was increased, the first 
eigen-value oscillated about a mean value whereas the higher 
order eigen-values converged. It was conjectured that this 
could possibly be due to 

(i) i yy ». I xx , leading to the ill-conditioning of the 

aystamjcf equations} and 

(11) the number of iterations in the NAG subroutine are 
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fixed but unknown and the user cannot exercise any 
control ox&r it* The eigen-values are the roots of 
the polynomial given by 

| A - * B i » 0 (2*7 *m) 

*•** 

The higher roots of the polynomial equation (2.7 *in). 
can be obtained by neglecting the lower degree terms 
while for finding out the lower eigen-values, in part- 
icular the lowest one, one has to consider all the 
terms* Consequently, higher roots ©an be computed 
with fewer number of iteratilns whereas the first few 
roots can be obtained only after a considerable number 
of iterations* Therefore, it is quite apt that the 
value of the first eigen-value is not a converged 
value within the number of Iterations fixed by the 
subroutine itself. 

Since thf number of iterations is beyond the control 
of the user, to mitigate the problem of , oscillating first 
natural frequency, one solves the Mi syatsm as an inverse eigen- 
value problem so^that the lowest mode becomes the highest mode 
of the new system. Thus, in lieu of solving the system of 
equations (2.7*1), one solves 



where M i 


Mk U ** £ 

_J 

^(n-i+t ) 


(2*7*n) 

.*** 


**•* 

The number of eigen-values ± for the new system of equations 
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is equal to the order of K or W» The eigen-values are the 
inverse of the square of i* h natural frequency^* Hqnce, ... 
the first (lowest) eigen-value ^ of the actual system., becomes 
the last (highest) eigen-value of the Inverse system 
(jif m ^—), The converged value of ^ is the reciprocal of 
th. converged value ot U K ' “ 13 verified that thi* indeed 
yi elds.. converged eig.gn— values with the increase in the number 


of elements* The number of frequencies and the corresponding 
eigen-vectors obtained is equal to the order of the partitioned 
stiffness or mass matrix and can be found out by subtracting 


thj» number of the speci f i ed 4< degrees of freedom from the total 
number of degrees of freedom of the approximated model* 


2*3*2 DATA AND RESULTS * 

The results have been obtained for the blade used in 


[12] so that the results can be 
for the aerofoil section are* 


*»* 


X 11 

as 

3*49634 x 1CT 11 m 4 

*«*. 

*22 

m 

2.79291 x 1C* 9 m 4 

*** 

1 

m 

0*1524 m 

r x 

m 

1*9304 x 10* 4 m 

r 

y 

as 

1*1938 x 10~ 3 m 

*** 

A 

ax 

5 2 

5*8967 x 10 m 


compared with it* The data 


Geometric 

properties 
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E * 213.7388 GPa 

**+ 

C * 9.1485978 N/m 
<? * 7858.3702 kg/m 3 


Water! al 
properties 


To test whether the given cross-section can be 
treated as a beam or not, the Zukovsky transformation 22 was 
made use of. The Zukovsky transformation, which maps an aero- 
foil-section in -plane to a circle in z-plane as shown in the 
Figure 2.3, is given by 

2 

£ * t + 5- where ? « ? + i*) 

z 1 (2.8, a) 

and z * x + iy 


Refer to Appendix I for details* 

... ... The transformation gives the chord length b and the 
maximum thickness h as 

b * 4(l 22 /A)^ 2 (2.8. b) 

and 

h » 3(2I t1 /A) 1/2 (2.8.C) 

From equations (2.8.b) and (2.8. c) one gets 
g* ss ^( A/I ** 5.53606 

and 

| * 3 tl 22 /2I 11 )1/2 * 8*426461 
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As the length, is large in comparison with the chord 
length and the maximum thickness, ^for the data considered, 
the blade can be treated as a beam. 

Results i For the data given above, the computed results 

upto fifth mode are given in the Table 2*1 and are compared 
X£th those obtained in the reference [J2j by (a) analytical 
method and (b) Ritz-Galerkin method. The physical identifi- 
cation of the frequencies was done by drawing the correspon- 
ding mode shapes (they ar.§ not shown here). The results were 
obtained for fourteen elements of equal length* 


... Table 2.1 

1 b 

Comparison of results for •£ - 5.53606, j- » 8.426461 


Method 


1 — — — — 

Descriptioid Frequency in Hz 

of 1 i -t r r~" 

modes }Fi rst {Second {First {Fi rst {Third 

{bending {bending {bending {torsion {bending 

{in YY Jin YY {in XX { {in YY 

{I } II ? HI { IV {V 

i .i, jt ... j- ■> 


Analytical results 02] 

Results by Ritz- 
Galerkin method 02] 


96.9 

606.5 

841.2 

1072.9 

1699*0 

96*5 

604*9 

842.0 

1089*4 

1 694*9 

96*8 

606.3 

843# 3 

1078*5 

1697*8 


Present results 
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2.4 CASE 2: 


Having studied the vibrations of untwisted cantilever 
beam of aerofoil cross-section, we now focus our attention on 
the pre-twisted cantilever beams of rectangular cross-section. 
For a rectangular section, the centroid and the flexural 
centre are coincident. So, the torsional mode is uncoupled; 
but the flapwise bending and chordwise bending are coupled 
because of the pre- twist. If the pre- twist is zero, both the 
bending modes also get decoupled. In this case, it is assumed 
that the effect of fibre-bending is not negligible. 

If b and h are the breadth and the height of the beam 
of rectangular cross-section, the expressions for the constants 
C.j , and Cg in equations (2.1. a) and (2.1.c) are given as 
follows: 


(i) 


n 


\ Gbh 3 {1 


192 h 1 

R K " 

it n=1,3,5... n 


tanh 


rut b -i 
2 h * 


(2. 9. a) 

where is the St. Venant torsional stiffness [19]. 


(ii) 


C 3 



( 2.9.b) 


where Cg is the torsional stiffness due to pre-twist 

[ 10 ] . 


(iii) 


Eb 3 h 3 
2 " 1 44“ 


(2.9. c) 
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where C 2 is the constant for the fibre-bending f20] • 

In this case, the tensors A and C of the equation 
(2*3*b) are given by 


* 


i * ij + EI mu £2 and £ * dia 9 L? a * ?A, ^cf] <2*10#a) 
where 


and 



A t * d ^ a ^-^av* > ^av* *^2^ 

cos sin 0 

sin — j*— —cos 0 


{ 2* 10*b) 


L 0 0 

* 

As the fibre-bending is considered, the maximum deri- 
vative of 9 in the functional is of second order* Hence the 
principal derivatives for 9 are ^ and 9* The convergence 
criteria states that 9 should also be expressed by Hermittan 



polynomials* P i , s. Hence, in equation (2.4.a), one gets 

(U® ) T * L U 1 %z m v i S? 3 * 


and 


S * 

m 


w * * v i - ISi „ ^ 

2 ^z 


'1 &z 


»'- r O ^ 0O 

U m , T .fe t/ fa fft ..- rT ~,--y 1 

ri ^T“ v ft J 


2 2z 


2 3Z 

(2*1t*a) 


P, p 2 0 0 0 0 P 3 P 4 0 0 0 0 

0 0 P t P 2 0 0 0 0 P 3 P 4 0 0 

0 0 0 0 P 1 P 2 0 0 0 0 P 3 P 4 


(2*1t*b) 
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The expressions for element stiffness matrix and mass matrix 
are the same as those in the previews ease* With the incor- 
poration of fibre-bending, the order of the element stiffness 
and mass matrices becomes 12 x 12 and can be evaluated by 
simplifying Integra tions^are discussed in the Section 2*3* 
Applying the geometric boundary conditions, 

*U t 2>V, 2>®. 

at z * 0 t Uj # 0, ^"35 * 0, Vj *s 0, g-— 1 «s 0, 0^ *e0, 11111 * 0 

(2.11 »c) 

which physically mean that the deflection and the deflection 
derivatives are zero at the fixed .end of the cantilever beam, 
one obtains the eigen-value problem on the lines similar to 
the Section 2*3. 

i> * , 

2.4.1 DATA A NO RESULTS * 

The results have been obtained for the blades used in 
riC , r 8J and Hit] so that the results can be compared with 
these cases. 


2.4.1. a Blade 1* The effect of the fibre-bending 1 


Data 1 b * 0.0254 m 


h * 0.001727 m 

*•* Geometric properties 

1 * 0.1524 m 

r« J rad 
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E = 206.85 GPa 

G = 82.74 GPa Material properties 

p = 7858.3702 kg/m 3 

•j 

As ^ = 6 and ^ = 14.706, the length is larger than the other 
dimensions and hence the blade can be treated as a beam. 

Results: For this data, the frequencies obtained are given in 

the Table 2.2 and their modes are identified as discussed in 
the previous case. The beam was discretized into 14 elements. 


Table 2.2 

Frequencies for ^ = o, |* = 14.706 and Y = ^ rad 


Mode 

No. 


Frequency in Hz 


— — — r — ■ r 

Bending in Y-Y ! Bending in X-X . 
direction (coupled); direction (couplea), 


Torsion (uncoupled) 


Present {Results of {Present 
results {ref. [18J {results 


{Results of {Present 
{ref. [18J {results 


{Results of 
{ref. L18J 


1 62.027 62.17 

2 310.005 310.18 

3 1232.751 1232.54 

4 2223.903 2223.51 


977.620 977.62 756.421 758.15 

5930.315 5926.49 2298.804 2349.61 

3928.691 4162.77 

5702.135 6242.28 


It is observed that the frequencies of bending-bending 
modes tally with the results given in the reference [l8], but 
the torsional frequencies do not agree with the results. This 
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is because in the reference [18], 

(i) the section has been aporoximated to the elliptic 

form and the corresoonding formula for the St. Venant 
torsional stiffness has been used and 
(ii ) the contribution of the pre-twist to the torsional 

stiffness is neglected and also it is assumed that the 
fibre-bending effect is negligible. 

2.4.1.b Blade 2: The effect of the pre- twist: 

Let us now study the effect of pre- twist on the vibra- 
tions of the beam with rectangular cross-section. For this, 
the data used in the reference [8] is considered. Using this 
data, the frequencies are found out for different values of 
pre-twist. The graph of reference [8] is shown in Figure 2. 4. a. 

Data: b = 0.0254 m 

h = 0.00635 m Geometric properties 
1 = 0.3048 m 

E = 206.85 GPa 

G = 82.74 GPa Material properties 

3 

P = 7858.7302 kg/m 

For this data, ^ = 12 and | = 4 and hence the -blade is treated 
as a beam. The frequencies for various values of pre-twist 
are tabulated with their modes identified and their trend 
with respect to the pre-twist in the Table 2.3. The beam is 
discretised into fifteen elements. 



Table 2.3 

The effect of f> re- twist on the frequencies for 
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It is observed that each m 0 de has got a specific 
trend of variation with respect to the angle-of-pre- twist. 

The first frequency corresponding to the first Y-Y bending 
mode increases with pre- twist. The second mode corresponding 
to the first X-X bending mode decreases with pre-twist and 
the fifth mode corresponding to the first torsional mode 
increases with the pre-twist. This is because the beam 
becomes stiffer for the Y-Y bending and the torsion whereas 
it becomes more flexible for the X-X bending at their first 
modes. The graphs in the Figure 2.4 depict the trends of the 
first nine freouencies with respect to the angle of pre- twist. 
These trends tally with those illustrated graphically in the 
reference [8], Figure 2. 4. a. 

2.4.1.C Blade 3: The effect of the pre- twist and the 

dimensions : 

For this blade, the data from the reference [17] is 
utilized. 

Data: h = 0.025 m Thickness 

E = 200.00 GPa 

= 0.3 (Poisson’s ratio) Material properties 
P = 7800 kg/m 3 

For this data and for four different geometries, the non- 
dimensionalized frequencies given by 
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where 


>1 ^Ph/D 


12(1 - |i ) 


( 2 . 12 ) 


are computed for the different angles of pre-twist and are 

tabulated in the Tables 2.4, 2.5., 2.6 and 2.7. The beams are 
discretised into seven elements. 

Table 2.4 

1 - i b - R 

F 1 * h = 5 


Pre- 

twist 

Y 

rad 


Y-Y 

ben- 

ding 


Tor- 

sion 


Non-dim ensionali zed frequency gl> 0 

hi ! iv ; v ! vi i vii i viii 

1 i i i 1- 

X-X ! Tor- ; Y-Y ! Tor- J Tor- \ Y-Y 

ben- ! sion { ben- J sion } sion j ben- 
ding ; J ding ; * | ding 


0 3.354 5.940 16.771 18.739 21.044 34.040 49.225 59.311 

it/12 3.36 3 6.154 1 6.087 18.615 21.98 9 31.5 34 45.221 58 . 487 

7i/o 3.37 1 6.410 1 4.743 1 9.38 9 24.006 32.8 46 47.10 3 56.908 

71 /4 3.379 6.677 13.432 20.196 26.440 34.213 49.062 54.722 

7t / 3 3.397 6.954 12.257 21.034 29.086 35.633 51 .099 52.315 


7t/2 3.448 7.539 1 0.36 3 22.8 0 4 34 . 38 3 38.6 3 2 55.399 48.041 
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Table 2.5 

|- 1 . | - 20 



0 3.354 5.364 19,693 21.030 34,889 52.348 59,067 67.082 


Vl2 3.466 9.453123,592 20.560 48,436 69.459 58.375 70*140 
V6 3.464 11,949 36.142 19.079 61,226 87.800 57.282 77.227 


*/4 3.395 14.167 42.850 17.211 72.590 104.097 56.123 86.339 


"V3 3.417 16.229 49.090 15.462 83.160 119,255 54.638 96.053 
V2 3*474 20.105 50.811 12,567 103,016 147.730 50.695 109,752 


Table 2.6 



Pre- 

twist 

r 

rad 



0 3.354 16,771 17,780 21,030 55.101 59,067 97,622 105*150 

tt/12 3*363 16,057 17.808 21.989 5 3.365 58,487 91,250 106.081 

v/6 3.371 14.743 17,898 24.006 54.138 56.908 91*712 108*357 

-"/4 3,379 1 3.430 17.996 26.440 54.432 54,722 92.210 110.079 

n /3 3.397 12.257 18.100 29.086 54,748 52.315 92*746 109*734 

n /2 3.449 10.363 18,333 34.383 55.453 48.041 93,939 1 06,366 
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Table 2.7 

* 3, | *» 20 


Pre- 

twist 

r 

rad 


Non-dim ensionalised frequency 


T 
t 
f 

'I 

1 Y-Y 
J ben- 
J ding 


T - ~T- — *T*~ — 

{ III } IV } V 

"I 1 i 

J Tor- 1 Tor- J Tor- 

{ sion { sion J sion 

i » » 

t t f 

Jt i — 


~r 

t 
* 

■\ — »- 
! Y-Y ! 

j ben— | 

} ding { 

j. l. 


II 


VI 


VII 


VIII 


-? 

J Tor- 
} sion 


Tor- 

sion 


Y-Y 

ben- 

ding 


0 

1.903 

3.354 

5.908 

10.466 

15.695 

8 

v— 

e 

O 

<N 

21.030 

59.067 

Vi 2 

7.269 

3.456 

21 .986 

37.246 

53.412 

70.821 

20,560 

58.375 

76 

10.308 

3*463 

31.179 

52.819 

75.744 

100.433 

19.079 

57.282 

"74 

12.813 

3*398 

38,756 

65*654 

94.151 

124.840 

17.211 

56.124 

73 

15.062 

3.420 

45,560 77.180 

110.680 

145.756 

15.462 

54.638 

72 

19*175 

3.467 

57.999 

98.253 

140.898 

186.824 

12.566 

50.595 


The trends of the frequencies versus the pre- twist for 
the four different aspects ratios are plotted in the graphs — * 
Figures 2.5, 2.6, 2*7 and 2.8 — corresponding to the Tables 
2*4, 2.5, 2.6 and 2.7 respectively* As surmised, some modes 
dither exhibit the trends contradictory to those in the refer- 
ence j^lTjor do not exhibit a definite trend, for the short blade 
is treated as a beam here. 











x,y,z, = Reference axesi1,2 = principal axes 
OO-i = Flexural axis ;PQ sCentroidal axis 
)f= Total angle of pretwist. 

FIG.2.1.a PRE-TWISTED CANTILEVER BEAM 


cfi centre of flexure; cm= centre of mass 
FIG2.1.b CROSS-SECTION OF THE BEAM 
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N B- The encircled numbers represent 
the element numbers. 

FIG. 2.2a DISCRETIZATION OF THE CANTILEVER 
BEAM INTO A FINITE NUMBER OF TWO-NODED 
ELEMENTS. 



eth element 

Z -Global co-ordinate, I -Local co-ordinate 
a -Half the length of the element 
Numbers inscribed in the triangle represent 
local node numbers. 

FIG. 2.2. b A TYPICAL TWO-NODED ELEMENT 






mi -twist a**o£ 


FIG. 2.4 


.a GRAPH OF 

FREQUENCY VERSUS -PRE-TWIST 
FOR L/B = 12 AND B/H = 4 
REFERENCE 8 




FIG. 



PRE-TWIST 
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CHAPTER 3 


VIBRATIONS OF BLADE AS A CANTILEVER PLATE 

3.1 BtAOE AS A PLATE ! 

When the length of the blade is short and Is compa- 
rable to Its cross-sectional dimensions, it tends to behave 
like a plate rather than a beam* If the blades are pre- 
twisted, they have geometric curvatures and torsion and hence 
behave like shells* When their pre-twist is zero, their 
geometric curvatures and torsion are zero and can be treated 
as a plate* A plate is defined as an object bounded by two 
flat surfaces* If the distance between these surfaces is 
very small as compared to the other dimensions of plate, it 
Is called a thin plate, else it is said to be a thick plate* 
For the plates, where the thickness is very small as compared 
to either the overall dimension or to the wavelength corres- 
ponding to the highest frequency of interest# the effects 
due to shear deformation and rotary inertia can be neglected 
without affecting the accuracy of the solution QC* 

This chapter deals with the study of free vibrations 
of a short blade treated as a uniform thin cantilever plate 
of rectangular cross-section, using FEW and neglecting rotary 
inertia and shear deformation* 

3.2 DERIVATION t 

For the thin plate shown in Figure 3*1, let w(x, y) 
be the displacement of a point on the middle surface along 
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x»axis* If we assume that the normal to the middle surface 
continues to be normal to it even after the deformation, 
then the displacements u and v of any arbitrary point 
(x, y, z) are given by 


u * z\7w i.e. u = z 


V = z 


3w 

Bx 

Bw 

By 


(3*1 • a ) 


The most significant strain and stress components are 


2 

£ *® zV V w i • e* £ = z ^ -a 

XX 

2 

*= z — o 
YY 9 y 2 

r * z JLSU 
b xy z BxB y 


(3.1. b) 


and 


H 1 1 1 _ ^ 


»i+ (1 -A) z V V W U 

— a2 _ 2 

w . /j d W) 


1 - M 

i.e. s'.. = 7^75 ^ ^ 

.2 ,/> 


XX 


^ „ Ez jfsu 

° xy 73^ 3x By 


(3.1.c) 


The stress-resultant over the thickness are defined 


by the relations 
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h/2 

M x * S z dz 

-h/2 xx 


h/2 

M y * $ z dz 

y -h/2 


yy 


h/2 

M xv “ * Jj’yv z dz 

xy -h/2 xy 

h/2 

Qy ** 5 tfl- dz 

x -h/2 Xz 


Q ss 

y 


h/2 

f dz 

-h/2 yz 


( 3.2*a) 


The stress- resultant M represents the bending moment per 
unit length on x-face about y-axis. Similarly m also repre- 
sents the bending moment per unit length but on y-face about 

x-£«fce* The quantity M denotes the twisting moment per 

xy 

unit length about x-axis on x-face (or about y-axis on 
y-faced* The shear forces per unit length along z-axis on 
x and y-faces are represented by Q x and Q y * The relation 
between the moment tensor 


M 

•m 


SB 


M M 

x xy 


M xy M y 


and the shear force vector 


a 


r -> 

Q 

M x 


iss 


1 
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is given by 


The expressions for jL, and C£ in terras of w ares 



u 

■ra 

*= D[/iv^ w + (l - /U) w w) 

(3.2.b) 

and 

a 

* D V(v 2 

w) 

(3.2.c) 

where 

n 

Eh 3 

7T" 

(3.2*d) 



12(1 -f/ 



D is called the flexural rigidity of the plate. 

The strain energy £l6j of the plate is given by 


V « |Jf(l -f()(V7w).(7Vw)+ / U(v 2 w) 2 JdA 
A 

n /v r ^2 2 -\2 2 -\2 *~\2 

*= o IT C ( — $) + ( — 5) + 2yU(r~^)(-~^) + 

0X By <?x d y 

2(1 J dx dv - (3 - 3) 

The expression for the kinetic energy [26]] of the plate is 
given as 

T = i $ h // (||*) dx dy (3.4) 

Applying Hamilton* s principle, 

*2 

<§ / (T - V) dt = 0 ( 3.5.a ) 

t* 


one obtains the governing equation of motion of the plate as 
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D \7 W + 


S h |^ 


{ 3*5»b) 


At the fixed end of the cantilever plate, the displacement 
and the displacement gradient are zero* Hence, the geome- 
tric boundary conditions at x * 0 are 

1 ) w m 0 and I~* *= 0 

o x 

called the essential boundary conditions and (3.5.c) 

ii) ~ « 0 and -SLJaL. * o 

11 ' 3y 'b x a y 


called the implied boundary conditions. 

At the free boundary, the bending moment, the twisting moment 
and the shear force, are all zero. That is, n = 0 and 
q . n « 0, where n is the unit outward normal to the boundary. 
Hence, the natural boundary conditions for a rectangular 
plate shown in Figure 3*t are 


i) 


At x = 1, 


A , u dfw 

_ +J U ™ 
3X dy 


d 2 W- *. 0 
BXB y 


0 


'*!» 

3vc 


c)\ 

- 

Sx S Y 


0 


b 


il ) At y w t 5 * 


3 2 w u 3 w 

"5 a*- r 2 

5y X 
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3 2 w 

" ■ » »» - ■ ■ i aes 

3x 3 y 


0 


3 3 w 3 3 w 

3x^0y 3y J 


0 


(3*5.d) 


At the governing equation is of the fourth order, 
two boundary conditions are needed at each boundary, whereas, 
in equation (3.5. d) three boundary conditions are obtained at 
each free boundary* In order to make the boundary conditions 
compatible, the boundary conditions for the shear force and 
the twisting moment are clubbed together £16”] and hence, the 
equivalent natural boundary conditions at the free edges of 
the rectangular plate are 


i) At x = 1, 


>2 >2 

+ . o 

3x z 3y 


o3 

— w + ( 2 — fJ . ) 
e>X 


3 3 W 
3x 9y S 


o 


ii) At y * ± | » 




(2 -jU) 


3 3 w 
3y 0x^ 


*5 


0 


C 3*5# ©) 
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3 * 3 FINITE ELEMENT FORMULATION : 

The closed form solution for a rectangular plate can 
be obtained explicitly if and only if the edges of the plate 
are simply supported. However, in the case of cantilever 
plate, closed form solution is difficult to obtain. There- 
fore, to study the vibrations of a cantilever plate, one has 
to make use of some numerical technique for the given plate(s) 
and here, the FEM has been opted for. 


Consider the expression for the strain energy in 
equation (3.3). The highest derivatives of w in the expres- 

A ^ 


o 2 

a w aw 


sion of the strain energy are — ^ and § ' x "jjf T * There 

3x 3y aw aw 


fore, the principal derivatives are w, and —• The 
compatibility condition states that these principal deriva- 
tives should be continuous across the inter-element boundaries. 
The highest derivative of w in boundary condition is of third 
order. Completeness condition states that the shape func- 
tions should be such that, as the element size reduces, the 
function and its derivatives upto the highest order accuracy 
in the variation of the strain energy are of finite values. 

This leads to the conclusion that the shape functions should 
be the product of at least cubic polynomials in t and T). 

It is noted that a mesh of rectangular elements is used 
to suit the convenience of the geometry and because of the ease 
of finding the compatible shape functions. Consider a four 
noded rectangular plate element as shown in Figure 3.2. If 
(x v ) is the global mid-point of the element, p and q 


, F18 S 35 
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are half the length and breadth of the element respectively, 
then a point within the element represented by the global 
co-ordinates (x, y) can be expressed in terms of the local 
co-ordinates as 



where -1 ^ ^ < 1 

and -1 < h ^ 1 

The displacement of any point within the element in 
terms of the nodal degrees of freedom and the shape functions 
is 


w = S . W e e lCOt 


(3.7) 


where W®, the nodal displacement vector, is given by 


(W e ) T 


3« 1 b 2 w 1 ti ^w 2 3W 2 ^ 2 
Sx dy 3x3 y W 2 'Six oY dx dy 


3 «l aw. . 3 w. A. 

nr - r[ . m y *3 v ur _ 4. M Sl 

3 "ax 3y Qx a y 4 3x ay 'dx'ciy 


( 3*8 • a ) 


Thus, there are sixteen degrees of freedom per element. The 
shape function vector is, 

S = 0 1 0 2 0 3 0 4 0 5 0 6 0 ? 0 Q 0 9 0 1Q 0 0 12 ^ 13 0 i4 0 15 0j 6 

(3.8*b) 

The expressions for 0^’s which satisfy the convergence crit- 


eria are: 
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For 1 = 1,2, 3, 4, 

i 

*4i-3 “ 16 (2 + (2 + 3 7 ± ^ - ^yj 3 ) 

^4i~2 * if ( ~?i -1 +iiX 2 +^ 3 )(2 + 37 i 9- 7 i ^ 3 ) 

**4i-1 * T§ (2 + 3^? - S^S 3 ) (- 7 ± - ? + + ^ 3 ) 

^41 = - i + i ± j 2 + ^ 3 ) (- y ± - yj + 7 ± ^ 2 +^ 3 ) 

(3.9) 

Substituting w = S . W e in equation (3.5), 

one gets the residue S e 


£ e = n(v 4 S) . W e . §h » 2 s . w e 

wkicl^ cab loe vc Written as 

f . Jr ^ „ ?4 s 
= ^ 5 jr + 


x 7* + 

3 x 3 y 2 3y' 


4 

9 S 7 e 
— ?5 . W + 2( t - 


- ?hg > 2 s . W® 


^ 4 S 

• l e 

d x d y 


( 3. 10*a) 


Applying Galerkin's method, 


Vi Vi 

J / S £ dx dy = 0 

x„ y 
e *e 


*W1 Vi r r J 4 S 

f * l D -l^ T ^ r~s 


2 4 S 9 4 S 

^ + 2M 

^x^ ^x^y^ 3 y 


9 4 s 


5 + ^ + 2(1 


- S h o) 2 SSdx dy)V = 0 


(3.10.b) 


Integrating by parts, the partial differential equation 
reduces to 
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K® ?J® - * «*• W e « R e 


(3*10. c) 


® SS ^ 2S are element stiffness matrix;, element mass 

matrix and element stress-resultant vector and are given by 


K® « D f " f' \ (1 -fJi ){W s).(V V s) + 
x e Ye 


*e+1 ^e+1. 


/( (v 2 s)(v 2 S) \ dx dy 


rrr ^ 2 S B 2 S } 2 S ^ 2 S tfs 


3 2 s Q 2 s 1 

+ 2(1 ~ M ) §Tay aTay j dx dy 


(3.l0.d) 


M 


@ 


3 h {J S S dx dy 


( 3. 10. e) 


R 


* *e+1\3 *^3 ^S\2 ^2 

W D I f s(— g + S . W ) — “( 2 — T g)(s . W ) 


'e 


3 x 3 x d y' 




^ p ~j e+1 

(1 -P'lfy&n ( 5 -S e >] dy 


x 


e+1 


+ D r T S&-* + 2 x)(S. W e ) - ^=(f-R +A ~~*)(S. W®) 

3 y 3 y c)x y 3y‘ t 3x 


e, . !i ( i 2 
y 

y. 


a 2 


25^2 _ / e+1 

- (1 -^JjxfcTy ( ***8 V dx 


(S.KUf) 


The element stiffness matrix and element mass-matrix 
are computed invoking Gauss-quadrature formula. 


Is 
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#* /• / \ **yA J 

{ S fix, y) dx dy * x: w.w. f(x 4 , y . ) (3*11) 

1 -1 i=1 j=l 1 J 1 J 


in which w^s are weight functions corresponding to the Gauss 
points (x^t Yj), ngx and ngy being numbers of Gauss points 
along x-direction and y-direction respectively* In order to 
have the limits of integration -1 to t as indicated in equa- 
tion (3*11 ) , one transforms the variables in (3*10. a) and 
(3*10#b) info natural co-ordinates using the transformation 
given by equation (3*6). From this transformation, one obtains 

dx dy = pq d^d^, 

1 3 2 

a? = 

(3.12) 

l 

and — - c i- 

and ^50 y pq 

Therefore, the element stiffness matrix and element mass matrix 
are obtained as 






JL “ 2(5 i 

Pfl 3^2 5g 2 


q ^ 


ih 

If} 2 


2il=El 


7>% 


c )% 


pq yj tedy 


■} 


d|d/? ( 3* 12*a) 


and 
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1 

^ h pq S 
-1 




( 3*t2*b) 


As the limits of integrations in the equations 
(3*12*a) and (3.12*b) are in accordance with those in the 
equation (3.11 ), the integrations can be evaluated by Gauss* 
quadrature formula. The element stiffness matrix and mass 
matrix, are symmetric positive definite and are of the order 
16 x 16. The element matrices of all the elements are assem- 
bled according to the standard procedure leading to the global 
stiffness and mass matrices which are partitioned after speci- 
fying the essential and the implied geometric boundary condi- 
tions from equation (3. 5. a) so as to obtain an eigen-value 
problem 

f K - 'X Ml W * 0 ( 3.1 3*a) 

L ss *= J — — 


To solve this system, i.e. to obtain eigen-values and 
eigen-vectors, the NAG subroutine titled F02AEF is used. Its 
details have been discussed in Chapter 2. As in the case of 
the beam, here also the inverse eigen-value problem is solved 
considering the system of equations to be 

(3.1 3*b) 


j= 0 where fA 


3*4 DATA AND RESULTS i 

The results have been obtained for the blade consi- 
dered in reference fl^J • The data is as follows* 
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h = 0.025 m 

P = 0.3 

H = 200.0 GPa 

P = 7800.0 kg/m 3 

For this cata, the non-dimensionali sed freouencies given by 

o' ph \ ^ / ^ 

“o = (3.14) 

are computed for four different geometries. The first eight 
non-dimensionali sed freouencies corresponding to each of the 
four geometries are given in the Table 3.1 below. 

Table 3.1 

Freouencies of the plate 
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of nodes. The results are fairly close to those given in the 
reference [17] but do not exactly tally with it because the 
convergence of the frequencies with the increase in the 
number of elements is comparatively slow in this case and 
hence it requires more number of elements to acquire the 
desired accuracy. The increase in the number of elements is 
not possible owing to the limitations of the computer memory 
space required to store five square matrices and a few arrays. 
Needless to say, more accurate results cannot be obtained. 


/ 



ABCD= Neutral plane 
Q = (x.y^o) Psfx.y^) 
w= displacement of Q 
u.v = displacements of P 
along x and Ydirections. 


FIG. 3.1 FLAT PLATE 
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n x (n x+ l) 


n x = number of elements in the x direction 
ny=. number of elements in the y direction 
N.B. The encircled numbers represent the element numbers 

FIG- 3-2. a DISCRETIZATION OF A CANTILEVER PLATE INTO 
A FINITE NUMBER OF FOUR-NODED ELEMENTS 

T) 

L p 


A (-u) A(vo 



p,q- half the element length and breadth, respectively. 
N.B. The numbers inscribed in the triangle respresent local 
node numbers 

F IG-3-2.b A TYPICAL FOUR-NODED PLATE ELEMENT 
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CHAPTER 4 
CONCLUSION 


4.1 BLADE AS A BEAM t 

Table 2*1 shows the natural frequencies corresponding 
to the coupled bending-bending- torsion modes of a long *= 
5*53606), untwisted, cantilever beam of asymmetric aerofoil 
cross-section. The results are in agreement with analytical 

results of the reference [12]* 

Table 2.2 shows the natural frequencies corresponding 
to the coupled bending-bending modes and uncoupled torsional 
modes of a long (§ * 12), pre- twisted and cantilever beam of 
the rectangular cross-section* The results are in agreement 
with those in the reference [18] • The frequencies correspon- 
ding to flapwise bending and torsion increase wfetb the 
increase in the angle-of-pre-tvdst while those corresponding 
to chordwise bending decrease* This is consistent with our 
expectations, as the area moment of inertia about X-X and 
the product of inertia of the section increase, the stiffness 
due to pre-twist and the effect of fibre-bending increase 
while the area moment of inertia about Y-Y decreases* ^ 

Table 2*3 also shows the results for a long beam ( B - 
12) and there is a definite trend as far as the variation of 

the frequency with respect to the pre-twist is 

On the other hand, the results of Table 2*4 are for 
a , h „rt baa* (fc - 1 « 3). Hera so»e torsional or flap*!.. 



bendina frequencies decrease, chordwise bending frequencies 
increase while some frequencies do not exhibit any definite 
trend. This leads to the conclusion that 
(i) Long blaoes can be analysed sufficiently accurately 
using a beam theory but not the short blades. 

(ii ) ror long blades there seems to be a definite trend 
for the variation of frequency with respect to pre- 
twist. 

4.2 BLADE AS A PLATF. ; 

Table 3.1 shows the natural frequencies of a short 

blade for various length to breadth and breadth to thicknes 

ratios. These frequencies are fairly close to those given 

in the reference [17] . It is observed that the frequencies 
2 

for the same j- ratio are the same irrespective of the ^ 
ratio. One can, therefore, draw the inference that 
(i) Short flat blades can be analysed satisfactorily 
using the plate theory, and 

(ii) The frequencies of short blades are independent of 
^ ratio as long as it is sufficiently small. 
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APPENDIX I 


THE ZUKOVSKY TRANSFORMATION 22 


X ** z + • a - 


The Zukovsky transformation 22 
2 

whore z «= x + iy and £ * % 4 iy (I*t) 


is a conformal mapping that transforms a circle in the z-plane 
into an aerofoil section in the £ -plane as shown in the 
Figure 2.3. 

For any point p on the circle in the z plane, 
r *» q e cos© + s sin© 4 t cosj) (1.2) 

For small values of e and , 


(1.3) 


(- 4 5) 2 and 

'q r 

(1 ^ |) 29(1 4 cos©) 4 2/3 sin© 

Substituting z «* r e i@ in the Zukovsky transformation, one 
gets 

X * 2q cos© and ^ ^ 

2 

y » 2q e( 1 4 cos©) sin© 4 2q /& sin © 

Using (1.4), the thickness, the area, the area moments 
of inertia and the product of inertia can be obtained as* 

h • 4q e sin©( 1 4 cos©) 



A 


4n e 




-71 


q e( 3e^ + 5 p 2 ) 




■en 


„ 4 

4 it q e 


2ti q e P 


(1.5) 


when p is very small compared to e, the inclination of the 
principal axes 1 and 2 with respect to the £ - r, axes is 
quite small, usually of the order of 2°. Then 


I 11 

CM 

Lim 
p -* 0 

ii 

JJLft 

M 

, 4 3 

on q e 

I 22 

rsJ 

Lim 

P - 0 

I r i r ) - 

, 4 

47iq e 

For 

this 

case 

, the 

expression for q and e are 


( 1 . 6 ) 


a 


(l M)^ 
K 22 / J 


(2I u /3I 22 ) 


1/2 


(1.7) 


The chord length of the blade comes out to be 
b = 4q : 


4(I 22 /A) 1/2 


( 1 . 8 ) 


The maximum thickness is obtained by maximising the 
expression hg of equation (1.5): 

h = 3V3 q e = 3(2I n /A) 1/2 (I,9) 
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appendix II 

householder 1 s decomposition 


Normally a non-singular symmetric matrix P is 
decomposed as 

LpJ * luJ T lq] [u] (ii, i) 

where [tj] i s orthogonal and [q] is positive definite. However, 
Instead of expressing this transformation as a product invol- 
ving two matrices, one can express it in terms of a matrix 
and a column vector by expressing [uJ as 

L'fl » [I]- 2 [w|{w( T (II. 2) 

where W is a normalised column vector. 

It can be found out by the relation 


* 0, 

W 2 - 1 4(1 +%*> ] 1/2 

rnm 





for k*=3, 4... n 


where 


s ** j ( 1 — $ ^ j ) (ll. 3) 

and the sign in the expression for W 2 and can be chosen to 
be that of P 12 * 
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Using ( IX* 1 ) and (11,2) the decomposition for P 

can be expressed as 

Lp] « [g] - 2}w} { wf [q] . 2[Q] {w] {wf 4 

4 { w i M T M \w] {wj T (11.4) 

Advantage® t 

i) For the higher order matrix P t less computer memory 
locations are required as in lieu of a matrix U of 

the order n x n # a column vector w of n components 
Is utilized, 

ii ) Computational economy and numerical accuracy can be 
achieved simultaneously* 


% 



